home *** CD-ROM | disk | FTP | other *** search
/ TeX 1995 July / TeX CD-ROM July 1995 (Disc 1)(Walnut Creek)(1995).ISO / graphics / gnuplot / binary.c < prev    next >
C/C++ Source or Header  |  1993-09-15  |  15KB  |  526 lines
  1. #ifndef lint
  2. static char *RCSid = "$Id: binary.c%v 3.50.1.16 1993/08/27 05:04:42 woo Exp $";
  3. #endif
  4.  
  5. /*
  6.  * The addition of gnubin and binary, along with a small patch
  7.  * to command.c, will permit gnuplot to plot binary files.
  8.  * gnubin  - contains the code that relies on gnuplot include files
  9.  *                     and other definitions
  10.  * binary      - contains those things that are independent of those 
  11.  *                     definitions and files
  12.  *
  13.  * With these routines, hidden line removal of your binary data is possible!
  14.  *
  15.  * Last update:  3/29/92 memory allocation bugs fixed. jvdwoude@hut.nl
  16.  *               3/09/92 spelling errors, general cleanup, use alloc with no
  17.  *                       nasty fatal errors
  18.  *               3/03/92 for Gnuplot 3.24.
  19.  * Created from code for written by RKC for gnuplot 2.0b.
  20.  *
  21.  * Copyright (c) 1991,1992 Robert K. Cunningham, MIT Lincoln Laboratory
  22.  *
  23.  */
  24. #include <stdio.h>
  25. #if !defined(apollo) && !defined(sequent) && !defined(u3b2) && !defined(alliant) &&!defined(sun386)
  26. #include <stdlib.h> /* realloc() */
  27. #else
  28. #include <sys/types.h> /* typedef long size_t; */
  29. extern char *realloc();
  30. #endif
  31. #include <math.h>
  32.  
  33. #if defined(MSDOS) && defined(__TURBOC__) && !defined(DOSX286)
  34. #include <alloc.h>        /* for farmalloc, farrealloc */
  35. #define SMALLMALLOC
  36. #endif
  37. #if defined(_Windows) && !defined(WIN32)
  38. #include <windows.h>
  39. #include <windowsx.h>
  40. #define farmalloc(s) GlobalAllocPtr(GHND,s)
  41. #define farrealloc(p,s) GlobalReAllocPtr(p,s,GHND)
  42. #define SMALLMALLOC
  43. #endif
  44. #ifdef sequent
  45. #include <sys/types.h>      /* unsigned long size_t; */
  46. #endif
  47.  
  48. #include "plot.h"   /* We have to get TRUE and FALSE */
  49.  
  50. float GPFAR *vector();
  51. float GPFAR *extend_vector();
  52. float GPFAR *retract_vector();
  53. float GPFAR * GPFAR *matrix();
  54. float GPFAR * GPFAR *extend_matrix();
  55. float GPFAR * GPFAR *retract_matrix();
  56. void free_matrix();
  57. void free_vector();
  58.  
  59. /* versions of alloc, realloc and free that work with segmented
  60.    architectures (yuk!) */
  61. char GPFAR *
  62. gpfaralloc(size, message)
  63.     unsigned long size;        /* # of bytes */
  64.     char *message;            /* description of what is being allocated */
  65. {
  66. #ifdef SMALLMALLOC
  67.     char GPFAR *p;                /* the new allocation */
  68.     char errbuf[100];        /* error message string */
  69.     p = farmalloc(size);
  70.     if (p == (char *)NULL) {
  71.         /* really out of memory */
  72.         if (message != NULL) {
  73.             (void) sprintf(errbuf, "out of memory for %s", message);
  74.             int_error(errbuf, NO_CARET);
  75.             /* NOTREACHED */
  76.         }
  77.         /* else we return NULL */
  78.     }
  79.     return(p);
  80. #else
  81.     return alloc(size, message);
  82. #endif
  83. }
  84.  
  85. char GPFAR *
  86. gpfarrealloc(p, size)
  87.     char GPFAR *p;            /* old pointer */
  88.     unsigned long size;        /* # of bytes */
  89. {
  90. #ifdef SMALLMALLOC
  91.     return farrealloc(p, size);
  92. #else
  93.     return realloc(p, (size_t)size);
  94. #endif
  95. }
  96.  
  97. void
  98. gpfarfree(p)
  99. char GPFAR *p;
  100. {
  101. #ifdef SMALLMALLOC
  102. #ifdef _Windows
  103. HGLOBAL hGlobal = GlobalHandle(SELECTOROF(p));
  104.     GlobalUnlock(hGlobal);
  105.     GlobalFree(hGlobal);
  106. #else
  107.     farfree(p);
  108. #endif
  109. #else
  110.     free(p);
  111. #endif
  112. }
  113.  
  114.  
  115. /* 
  116.  * This routine scans the first block of the file to see if the file is a 
  117.  * binary file.  A file is considered binary if 10% of the characters in it 
  118.  * are not in the ascii character set. (values < 128), or if a NUL is found.
  119.  * I hope this doesn't break when used on the bizzare PC's.
  120.  */
  121. int
  122.   is_binary_file(fp)
  123. register FILE *fp;
  124. {
  125.   register int i,len;
  126.   register int odd;                /* Contains a count of the odd characters */
  127.   long where;
  128.   register unsigned char *c;
  129.   unsigned char buffer[512];
  130.  
  131.   if((where = ftell(fp)) == -1){ /* Find out where we start */
  132.     fprintf(stderr,"Notice: Assuming unseekable data is not binary\n");
  133.     return(FALSE);
  134.   }
  135.   else {
  136.     rewind(fp);
  137.  
  138.     len = fread(buffer,sizeof(char),512,fp);
  139.     if (len <= 0)                      /* Empty file is declared ascii */
  140.       return(FALSE);
  141.  
  142.     c = buffer;
  143.  
  144.     /* now scan buffer to look for odd characters */
  145.     odd = 0;
  146.     for (i=0; i<len; i++,c++) {
  147.       if (!*c) {              /* NUL _never_ allowed in text */
  148.     odd += len;
  149.     break;
  150.       }
  151.       else if ((*c & 128) ||/* Meta-characters--we hope it's not formatting */
  152.            (*c == 127)|| /* DEL */
  153.            (*c < 32 && 
  154.         *c != '\n' && *c != '\r' && *c != '\b' &&
  155.         *c != '\t' && *c != '\f' && *c != 27 /*ESC*/))
  156.     odd++;
  157.     }
  158.   
  159.     fseek(fp,where,0); /* Go back to where we started */
  160.  
  161.     if (odd * 10 > len)             /* allow 10% of the characters to be odd */
  162.       return(TRUE);
  163.     else
  164.       return(FALSE);
  165.   }
  166. }
  167. /*========================= I/O Routines ================================
  168.   These may be useful for situations other than just gnuplot.  Note that I 
  169.   have included the reading _and_ the writing routines, so others can create 
  170.   the file as well as read the file.
  171. */
  172.  
  173. /*
  174.   This function reads a matrix from a stream
  175.  
  176.   This routine never returns anything other than vectors and arrays
  177.   that range from 0 to some number.  
  178.  
  179. */
  180. #define START_ROWS 100/* Each of these must be at least 1 */
  181. #define ADD_ROWS 50
  182. int
  183.   fread_matrix(fin,ret_matrix,nr,nc,row_title,column_title)
  184. FILE *fin;
  185. float GPFAR * GPFAR * GPFAR *ret_matrix,GPFAR * GPFAR * row_title, GPFAR * GPFAR *column_title;
  186. int *nr,*nc;
  187. {
  188.   float  GPFAR * GPFAR *m, GPFAR *rt, GPFAR *ct;
  189.   register int num_rows = START_ROWS;
  190.   register int num_cols;
  191.   register int current_row = 0;
  192.   register float  GPFAR * GPFAR *temp_array;
  193.   float fdummy;
  194.   
  195.   fread(&fdummy,sizeof(fdummy),1,fin);
  196.   num_cols = (int)fdummy;
  197.   
  198.   /* 
  199.     Choose a reasonable number of rows,
  200.     allocate space for it and continue until this space
  201.     runs out, then extend the matrix as necessary.
  202.     */
  203.   ct = vector(0,num_cols-1);
  204.   fread(ct,sizeof(*ct),num_cols,fin);
  205.  
  206.   rt = vector(0,num_rows-1);
  207.   m = matrix(0,num_rows-1,0,num_cols-1);
  208.  
  209.   while(fread(&rt[current_row], sizeof(rt[current_row]), 1, fin)==1){ 
  210.     /* We've got another row */
  211.     if(fread(m[current_row],sizeof(*(m[current_row])),num_cols,fin)!=num_cols)
  212.       return(FALSE);      /* Not a True matrix */
  213.  
  214.     current_row++;
  215.     if(current_row>=num_rows){ /* We've got to make a bigger rowsize */
  216.       temp_array = extend_matrix(m,0,num_rows-1,0,num_cols-1,
  217.                  num_rows+ADD_ROWS-1,num_cols-1);
  218.       rt = extend_vector(rt,0,num_rows-1,num_rows+ADD_ROWS-1);
  219.       
  220.       num_rows+= ADD_ROWS;
  221.       m = temp_array;
  222.     }
  223.   }
  224.   /*  finally we force the matrix to be the correct row size */
  225.   /*  bug fixed. procedure called with incorrect 6th argument. jvdwoude@hut.nl */
  226.   temp_array = retract_matrix(m,0,num_rows-1,0,num_cols-1,current_row-1,num_cols-1);
  227.   /* Now save the things that change */
  228.   *ret_matrix = temp_array;
  229.   *row_title = retract_vector(rt, 0, num_rows-1, current_row-1);
  230.   *column_title = ct;
  231.   *nr = current_row;/* Really the total number of rows */
  232.   *nc = num_cols;
  233.   return(TRUE);
  234. }
  235.  
  236. /* This writes a matrix to a stream 
  237.    Note that our ranges are inclusive ranges--and we can specify subsets.
  238.    This behaves similarly to the xrange and yrange operators in gnuplot
  239.    that we all are familiar with.
  240. */
  241. int
  242.   fwrite_matrix(fout,m,nrl,nrh,ncl,nch,row_title,column_title)
  243. register FILE *fout;
  244. register float  GPFAR * GPFAR *m, GPFAR *row_title, GPFAR *column_title;
  245. register int nrl,nrh,ncl,nch;
  246. {
  247.   register int j;
  248.   float length;
  249.   register int col_length;
  250.   register int status;
  251.   float  GPFAR *title = NULL;
  252.  
  253.   length = col_length = nch-ncl+1;
  254.  
  255.   if((status = fwrite((char*)&length,sizeof(float),1,fout))!=1){
  256.     fprintf(stderr,"fwrite 1 returned %d\n",status);
  257.     return(FALSE);
  258.   }
  259.   
  260.   if(!column_title){
  261.     column_title = title = vector(ncl,nch);
  262.     for(j=ncl; j<=nch; j++)
  263.       title[j] = j;
  264.   }
  265.   fwrite((char*)column_title,sizeof(float),col_length,fout);
  266.   if(title){
  267.     free_vector(title,ncl,nch);
  268.     title = NULL;
  269.   }
  270.  
  271.   if(!row_title){
  272.     row_title = title = vector(nrl,nrh);
  273.     for(j=nrl; j<=nrh; j++)
  274.       title[j] = j;
  275.   }
  276.     
  277.   for(j=nrl; j<=nrh; j++){
  278.     fwrite((char*)&row_title[j],sizeof(float),1,fout);
  279.     fwrite((char*)(m[j]+ncl),sizeof(float),col_length,fout);
  280.   }
  281.   if(title)
  282.     free_vector(title,nrl,nrh);
  283.  
  284.   return(TRUE);
  285. }
  286.  
  287. /*===================== Support routines ==============================*/
  288.  
  289. /******************************** VECTOR *******************************
  290.  *       The following routines interact with vectors.
  291.  *
  292.  *   If there is an error we don't really return - int_error breaks us out.
  293.  *
  294.  *   This subroutine based on a subroutine listed in "Numerical Recipies in C",
  295.  *   by Press, Flannery, Teukoilsky and Vetterling (1988).
  296.  *
  297.  */
  298. float GPFAR *vector(nl,nh)
  299.      register int nl,nh;
  300. {
  301.   register float GPFAR *vec;
  302.  
  303.   if (!(vec = (float GPFAR *)gpfaralloc((unsigned long) (nh-nl+1)*sizeof(float),NULL))){
  304.     int_error("not enough memory to create vector",NO_CARET);
  305.     return NULL;/* Not reached */
  306.   }
  307.   return (vec-nl);
  308. }
  309. /* 
  310.  *  Free a vector allocated above
  311.  *
  312.  *   This subroutine based on a subroutine listed in "Numerical Recipies in C",
  313.  *   by Press, Flannery, Teukoilsky and Vetterling (1988).
  314.  *
  315.  */
  316. void 
  317.   free_vector(vec,nl,nh)
  318. float  GPFAR *vec;
  319. int nl,nh;
  320. {
  321.   gpfarfree((char GPFAR *)(vec+nl));
  322. }
  323. /************ Routines to modify the length of a vector ****************/  
  324. float  GPFAR *
  325.   extend_vector(vec,old_nl,old_nh,new_nh)
  326. float  GPFAR *vec;
  327. register int old_nl,old_nh,new_nh;
  328. {
  329.   register float  GPFAR *new_v;
  330.   if(!(new_v = (float GPFAR *)gpfarrealloc((void*)(vec+old_nl),
  331.                        (unsigned long)(new_nh-old_nl+1)*sizeof(float)) )){
  332.     int_error("not enough memory to extend vector",NO_CARET);
  333.     return NULL;
  334.   } 
  335.   return new_v-old_nl;
  336. }
  337.  
  338. float  GPFAR *
  339.   retract_vector(v,old_nl,old_nh,new_nh)
  340. float  GPFAR *v;
  341. register int old_nl,old_nh,new_nh;
  342. {
  343.   register float GPFAR *new_v;
  344.   if(!(new_v = (float GPFAR *)gpfarrealloc((void*)(v+old_nl),
  345.                                (unsigned long)(new_nh-old_nl+1)*sizeof(float)))){
  346.     int_error("not enough memory to retract vector",NO_CARET);
  347.     return NULL;
  348.   }
  349.   return new_v-old_nl;
  350. }
  351. /***************************** MATRIX ************************
  352.  *
  353.  *       The following routines work with matricies
  354.  *
  355.  *      I always get confused with this, so here I write it down:
  356.  *               for nrl<= nri <=nrh and
  357.  *               for ncl<= ncj <=nch
  358.  *  
  359.  *   This matrix is accessed as:
  360.  *   
  361.  *     matrix[nri][ncj];
  362.  *     where nri is the offset to the pointer to a vector where the
  363.  *     ncjth element lies.
  364.  * 
  365.  *   If there is an error we don't really return - int_error breaks us out.
  366.  *
  367.  *   This subroutine based on a subroutine listed in "Numerical Recipies in C",
  368.  *   by Press, Flannery, Teukoilsky and Vetterling (1988).
  369.  *
  370.  */
  371. float 
  372.    GPFAR * GPFAR *matrix(nrl,nrh,ncl,nch)
  373. register int nrl,nrh,ncl,nch;
  374. {
  375.   register int i;
  376.   register float GPFAR * GPFAR *m;
  377.  
  378.   if (!(m = (float GPFAR * GPFAR *)gpfaralloc((unsigned long)(nrh-nrl+1)*sizeof(float GPFAR *),NULL))){
  379.     int_error("not enough memory to create matrix",NO_CARET);
  380.     return NULL;
  381.   }
  382.   m -= nrl;
  383.  
  384.   for (i=nrl; i<=nrh; i++)
  385.     {
  386.       if (!(m[i] = (float GPFAR *) gpfaralloc((unsigned long)(nch-ncl+1)*sizeof(float),NULL))){
  387.     free_matrix(m,nrl,i-1,ncl,nch);
  388.     int_error("not enough memory to create matrix",NO_CARET);
  389.     return NULL;
  390.       }
  391.       m[i] -= ncl;
  392.     }
  393.   return m;
  394. }
  395. /* 
  396.  * Free a matrix allocated above
  397.  *
  398.  *
  399.  *   This subroutine based on a subroutine listed in "Numerical Recipies in C",
  400.  *   by Press, Flannery, Teukoilsky and Vetterling (1988).
  401.  *
  402.  */
  403. void 
  404.   free_matrix(m,nrl,nrh,ncl,nch)
  405. float  GPFAR * GPFAR *m;
  406. unsigned nrl,nrh,ncl,nch;
  407. {
  408.   register int i;
  409.  
  410.   for (i=nrl; i<=nrh; i++) 
  411.     gpfarfree((char GPFAR *) (m[i]+ncl));
  412.   gpfarfree((char GPFAR *) (m+nrl));
  413. }
  414. /*
  415.   This routine takes a sub matrix and extends the number of rows and 
  416.   columns for a new matrix
  417. */
  418. float GPFAR * GPFAR *extend_matrix(a,nrl,nrh,ncl,nch,srh,sch)
  419.      register float  GPFAR * GPFAR *a;
  420.      register int nrl,nrh,ncl,nch;
  421.      register int srh,sch;
  422. {
  423.   register int i;
  424.   register float GPFAR * GPFAR *m;
  425.  
  426.   /*  bug fixed. realloc() called with incorrect 2nd argument. jvdwoude@hut.nl */
  427.   if(!(m = (float GPFAR * GPFAR *)gpfarrealloc((void*)(a+nrl),(unsigned long)(srh-nrl+1)*sizeof(float GPFAR *)) )){
  428.     int_error("not enough memory to extend matrix",NO_CARET);
  429.     return NULL;
  430.   }
  431.  
  432.   m -= nrl;
  433.  
  434.   if(sch != nch){
  435.     for(i=nrl; i<=nrh; i++)
  436.       {/* Copy and extend rows */
  437.     if(!(m[i] = extend_vector(m[i],ncl,nch,sch))){
  438.       free_matrix(m,nrl,nrh,ncl,sch);
  439.       int_error("not enough memory to extend matrix",NO_CARET);
  440.       return NULL;
  441.     }
  442.       }
  443.   }
  444.   for(i=nrh+1; i<=srh; i++)
  445.     {
  446.       if(!(m[i] = (float GPFAR *) gpfaralloc((unsigned long) (nch-ncl+1)*sizeof(float),NULL))){
  447.     free_matrix(m,nrl,i-1,nrl,sch);
  448.     int_error("not enough memory to extend matrix",NO_CARET);
  449.     return NULL;
  450.       }
  451.       m[i] -= ncl;
  452.     }
  453.   return m;
  454. }
  455. /*
  456.   this routine carves a large matrix down to size
  457. */
  458. float GPFAR * GPFAR *retract_matrix(a,nrl,nrh,ncl,nch,srh,sch)
  459.      register float  GPFAR * GPFAR *a;
  460.      register int nrl,nrh,ncl,nch;
  461.      register int srh,sch;
  462. {
  463.   register int i;
  464.   register float  GPFAR * GPFAR *m;
  465.  
  466.   for(i=srh+1; i<=nrh; i++) {
  467.     free_vector(a[i],ncl,nch);
  468.   }
  469.  
  470.   /*  bug fixed. realloc() called with incorrect 2nd argument. jvdwoude@hut.nl */
  471.   if(!(m = (float GPFAR * GPFAR *)gpfarrealloc((void*)(a+nrl), (unsigned long)(srh-nrl+1)*sizeof(float GPFAR *)) )){
  472.     int_error("not enough memory to retract matrix",NO_CARET);
  473.     return NULL;
  474.   }
  475.  
  476.   m -= nrl;
  477.  
  478.   if(sch != nch){
  479.     for(i=nrl; i<=srh; i++)       
  480.     if(!(m[i] = retract_vector(m[i],ncl,nch,sch))){ {/* Shrink rows */
  481.       free_matrix(m,nrl,srh,ncl,sch);
  482.       int_error("not enough memory to retract matrix",NO_CARET);
  483.       return NULL;
  484.     }
  485.       }
  486.   }
  487.  
  488.   return m;
  489. }
  490.  
  491. float 
  492.    GPFAR * GPFAR *convert_matrix(a,nrl,nrh,ncl,nch)
  493. float GPFAR *a;
  494. register int nrl,nrh,ncl,nch;
  495.  
  496. /* allocate a float matrix m[nrl...nrh][ncl...nch] that points to the
  497. matrix declared in the standard C manner as a[nrow][ncol], where 
  498. nrow=nrh-nrl+1, ncol=nch-ncl+1.  The routine should be called with
  499. the address &a[0][0] as the first argument.  This routine does
  500. not free the memory used by the original array a but merely assigns
  501. pointers to the rows. */
  502.  
  503. {
  504.   register int i,j,ncol,nrow;
  505.   register float GPFAR * GPFAR *m;
  506.  
  507.   nrow=nrh-nrl+1;
  508.   ncol=nch-ncl+1;
  509.   if (!(m = (float GPFAR * GPFAR *)gpfaralloc((unsigned long)(nrh-nrl+1)*sizeof(float GPFAR *),NULL))){
  510.       int_error("allocation failure in convert_matrix()",NO_CARET);
  511.       return NULL;
  512.   }
  513.   m -= nrl;
  514.  
  515.   m[nrl]=a-ncl;
  516.   for(i=1,j=nrl+1;i<=nrow-1;i++,j++) m[j]=m[j-1]+ncol;
  517.   return m;
  518. }
  519.  
  520. void free_convert_matrix(b,nrl,nrh,ncl,nch)
  521. float GPFAR* GPFAR *b;
  522. register int nrl,nrh,ncl,nch;
  523. {
  524.     free((char*) (b+nrl));
  525. }
  526.